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Abstract 

A deterministic version of the Bak-Sneppen model is studied. The 
role of the Lyapunov spectrum in the onset of scale-free behavior is 
established and avalanches are interpreted as return times to a zero- 
measure set. The problem of accurate determination of the scaling 
exponents near the critical barrier is addressed using a characteristic 
function approach. A general equation for return times to a small mea- 
sure set is established and used to obtain information on the nature 
of Bak-Sneppen dynamics. 

1 Introduction 

In Ref. P several structure-generating mechanisms, in multi-agent dynamical 
systems, have been analyzed. A structure in a dynamical system is character- 
ized by the emergence of phenomena at scales different from those typical of 
the dynamics of the agents, when they evolve in isolation. A space- structure 
corresponds to phenomena operating at a scale much larger than the typical 
size of one agent and a time- structure is a phenomenon with a time scale 
much longer than the cycle time of the agents. 

An important role in the generation of time-structures is played by modi- 
fications to the Lyapunov spectrum, occurring either as a result of the agents' 
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interaction or of a change of parameters. A new structure is created each 
time one Lyapunov exponent crosses zero from above. In particular, in one 
of the mechanisms described in |T1, the joint effect of chaos in the individual 
agent evolution law and extremal dynamics, implies that all Lyapunov expo- 
nents tend to O"*" in the infinite agents limit. This leads to a characteristic 
behavior free of time scales. It has been conjectured in Ref. PP that it is this 
mechanism that is behind many of the dynamical manifestations of what has 
been called self-organized criticality (SOC). Here, this conjecture is pursued 
by analyzing a situation |2J where breakdown of self-organized criticality is 
exhibited. It is shown that the breakdown of SOC may in fact be understood 
from the behavior of the Lyapunov spectrum. 

There is some controversy concerning a rigorous definition of SOC. Here 
only the absence of natural time scales is emphasized, without much concern 
about space scales, nature of the driving, separation of time scales and other 
relevant issues useful for a precise characterization of SOC. 

Absence of time scales, as a property of dynamics, is naturally related 
to the Lyapunov spectrum. Time scales disappear whenever the Lyapunov 
exponents vanish. This leads in a natural way to power laws for both time 
and space correlation. However, scaling of the avalanches, that is, a power 
law for the return times to the self-organized state, is a subtler effect. This is 
easy to understand when the multi-agent system is formulated as a measure- 
preserving dynamical system. If the self-organized state A, which serves 
as reference for the counting of return times {k), is a non-zero measure set 
{fi (A) 7^ 0), the measure /i (A) itself serves as a natural time scale. Then, the 
large time behavior of the return times distribution would be dominated by an 
exponential factor exp {—ku (/i)). Therefore, for cases that fit in the ergodic 
dynamical systems setting, power laws may occur only if the reference set 
has vanishing measure. In the fi {A) —>■ limit the return times distribution 
is dominated by the pre-factor that multiplies the exponential and this one 
might be or not be a power law. 

When /i is an ergodic measure, by Kac's lemma, the mean return time to 
a set A is l/fi{A). Therefore when fi{A) = 0, the numerical evaluation of 
the return time (avalanche) law has to be carried out for a set slightly larger 
than A. This implies that it may be difficult to disentangle the prefactor 
dependence from the exponential one, leading to some uncertainty about the 
exact values of the numerically measured scaling exponents. In particular 
because, as seen in Sect. 2, the exponent factor u (/i) may be a non-trivial 
function of the measure. Instead of obtaining the return times distribution 
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from the histogram, a more robust way is described in Sect. 2, based on the 
construction of the characteristic function. 

Treating the muhi-agent system as a deterministic measure-preserving 
dynamical system, tools from ergodic theory may be applied. An important 
issue is the characterization of asymptotic sets and invariant measures as well 
as the nature of the relevant attractors and repellers of the scale-free system. 
These questions are discussed for a deterministic version of the Bak-Sneppen 
model jHj. The nature of the self-organized state as a zero-measure subset of 
the invariant measure is put into evidence. Then, a probabilistic avalanche 
equation, written for the return times, is used to extract the non-trivial 
nature of the avalanche process from simulation data. 

2 The deterministic Bak-Sneppen model 

The original Bak-Sneppen model jS] |3] may be converted into a deterministic 
dynamical system by defining 

Xi {t+l) = Ti [x^ X, (t) + (l - (x)) / (xi (t)) (1) 

where x = {xi} is the vector of agent coordinates and / a deterministic 
pseudo-random generator, for example 

/ (xj) = kx mod.l 

k = 2,3, ■ ■ ■ . Fj is a function which is nearly zero if i corresponds to 

the agent with the minimum Xj value or to one of its 2n„ neighbors and is 
nearly one otherwise. In Ref. PP the following function was proposed^ 

j=i+nv / \ 

rS'^fe)= n 1-11(1 + ^""^""^^^)" (2) 

j=i-nv \ /^i / 

which, for large a, satisfies the above conditions. 
Here, instead, the function 

r?' (j) - n 1 - (3) 

j=i—nY \ Z-^l=l / 

^Notice that in Ref.p there is a typing mistake in the definition of this fmiction 
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will be used, which has a similar behavior for T — > 0^. Using a similar 
function in a stochastic model Headj2] has shown the breakdown of scale- 
free behavior for finite non-zero T. By using t[^'^ (^x^ in the deterministic 

model one may compare the Lyapunov spectrum analysis with the numerical 
results of Head. In the T — limit both models are equivalent to the original 
Bak-Sneppen model. 
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Figure 1: Approximate Lyapunov spectrum for several T values in the de- 
terministic Bak-Sneppen model 

Fig.l shows an approximate calculation of the evolution of the Lyapunov 
spectrum for successively decreasing T values with agents. The results 
were obtained for N = 100, 1000 time steps and 2ny = 2. The spread in the 
T = limit, around the exact value for N agents (-^ log 2 = 0.02 in this case) , 
results from the finite number of time steps. The deterministic system de- 
fined by (III)-® closely corresponds to the stochastic system studied in Ref. [2]. 
Therefore the behavior of the Lyapunov spectrum explains why the distribu- 
tion of jump sizes and activation times obeys power laws only in the T — > 
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limit. In Ref. [2] the loss of criticality for finite T is interpreted in terms of 
the correlated or uncorrelated nature of the active site jumps. The Lyapunov 
spectrum interpretation is clearer. Notice also that it is only in the N oo 
limit that all the Lyapunov exponents reach O"''. It is only in this limit that 
all time scales disappear JJ. 

One of the most remarkable features of the Bak-Sneppen model is the ex- 
istence of a sharp probability distribution for the agents coordinates above a 
barrier b ~ 0.667. However, this sharp (full-measure) distribution is merely a 
fixed point for the one-agent marginal, which is the one- dimensional projec- 
tion of a zero-measure subset of the invariant measure in the N— dimensional 
space [T]. The fact that in the A^— dimensional space the self-organized set to 
which the system returns after each avalanche is of zero measure is consistent 
with Kac's lemma. Kac's lemma states that for an ergodic invariant measure 
fi the average return time to a set A of measure fi (A) is 1/fi {A). Therefore 
for a scaling law p (A;) ~ l/A;'^, r < 2 implies fi (A) = 0. 

To exhibit the zero-measure nature of the Bak-Sneppen state consider the 
distance process d defined by 



Fig. 2 shows the probability distribution of the distance process for sev- 
eral values of the barrier b. One sees that it is for b around 0.67 that the 
neighborhood of the d = point becomes of zero measure. 

The Bak-Sneppen self-organized state is an A^— dimensional hypercube of 
volume (1 — 0.667)^. This set has repelling directions corresponding to the 
agents that are active and neutral directions for all others. Not being an 
invariant set it falls outside the usual definition of "weak repeller". It has 
been called a "ghost weak repeller" in Ref. jl]. 

As mentioned before, the zero measure of this "repeller" makes the direct 
measurement of the distribution law of avalanches a delicate matter. On the 
one hand the barrier value defining the avalanches must be placed close to 
the critical barrier to avoid the exponential finite-measure effects. However 
because the average size of the avalanches is (A), the closer we are to the 
critical barrier the worse the statistics becomes. A bin size has to be chosen 
to fit the avalanches data to a power law and the accuracy of the exponent 
is found to be sensitive to the bin size. 

Logarithmic binning may improve the situation. However a more robust 
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Figure 2: Probability density for the distance process 



method is based on the construction of the characteristic function 
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from the data. For each avalanche size k, the probability density may then 
be obtained by numerical computation of the inverse Fourier transform 
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p(k) = — C (x) e-^'^^rfx 
27r 



(6) 



In this way all the data is used for each k instead of only the events in the 
neighborhood of k. 

For sparse data this is a very robust way to construct p (/c), which is no 
longer affected by the bin size choice. However in our case we are not only 
interested in an accurate determination of p{k) but also on the extraction of 
the power law prefactor that multiplies the measure- dependent exponential. 
In cases where the exponent of the exponential is also a non-trivial function 
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of the measure it is extremely difficult to separate the effect of the two terms. 
This is especially true if the measure-dependent exponent decreases with the 
measure (see the example below). 

More accurate statements, concerning the nature of the prefactors, are 
obtained by comparing the characteristic function constructed from the data 
with the characteristic functions for trial distributions 
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with c and z/q, (yu)obtained from normalization and Kac's lemma, (k) = ^ . 
This was done for several values of the measure, the results being displayed in 
Fig. 3 (for 2ny = 2). Both the exponential factor i/ and the scaling exponent 
a depend on the measure of the avalanche reference set. 



10" 



10 ■ 



10 



10" 

















Or " ' 











10 



10^ 10 



10 



1 

0.9 
0.8 

L 

r 

0.7 
0.6 
0.5 





G- ^ ^ 

^ ~ e - 





























10" 



10"" 10 



10"' 



Figure 3: The exponential argument i/(/x) and scaling exponent «(//) 

This puts into evidence the unavoidable uncertainties in the direct eval- 
uation of the scaling exponent, in particular because even very close to the 
critical barrier the scaling exponent still shows a noticeable dependence on 
the measure of the reference set. In the Appendix the characteristic functions 
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associated to the probability distributions in Eq.((7)) are discussed in detail 
as well as some other results concerning return times to small measure sets. 

The discussion above refers to the problem of direct determination of 
the scaling exponents. However, with the additional assumption, as done 
by several authors, of a scaling form for p{k) near the critical barrier, an 
estimate of the fi —>■ value may be obtained. This however relies on the 
validity of the scaling assumption. Assuming that close to yU = 

p{k,^,) = k-'^f{k'^^) (8) 

one obtains 



{k) ~ /i 



a-2 

s 

a-3 



2\ a — i (^) 



Then, from Kac's lemma, one obtains 



a 

a-3 



and from the numerical data for = 2 (Fig. 4a), ~ 2.07, leading to 
a ^ 1.067, s ~ 0.93. 

Another exponent relates the measure to the position of the barrier 

with (Fig. 4b) r] ^ 2.55. 



3 The avalanche equation 

Assuming a differentiable dependence of the avalanches on the position of 
the barrier, Maslov derived an equation for the probability density of 
avalanche sizes in the Bak-Sneppen model. Here, a similar equation will be 
obtained for the return times to a measure fi set. This equation is then ap- 
plied to the Bak-Sneppen model. Particular attention is paid to the avalanche 
merging kernel, which in the Maslov equation is usually assumed to scale like 
a power in the critical region, but that may actually have a much more com- 
plex form. 

Assume that a set parametrization may be chosen in such a way that set 
measures change in a differentiable manner. Differentiating with respect to 
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Figure 4: Data points used to infer the exponents in (A;^) ~ /i » and n ~ 

the measure, one obtains the following equation 

27 {k, lJ^)dfx J - (j) 7 U, df^P {k ~ j \ j) (10) 

where 

P(T^^^ = A;,/i(A)) =pik) 

Eq. ljlOj) is derived by shrinking the measure of the set A by the amount 
and accounting for the number of avalanches of size k that become avalanches 
of larger size and for the number of mergers of smaller avalanches leading to 
new A;— avalanches. 7 {k, /i) d^ is the probability that the first or the last step 
of a fc— avalanche iteration falls in dfi. The first term in Eq.()10|) accounts 
both for a fc— avalanche becoming a larger avalanche and for the change in 

the total number of events from to ^1 + . The second term accounts 

for the probability for a j— avalanche and a (/c — j) —avalanche to merge into 
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a A;— avalanche. p{k — j \ j) is the conditional probability to have a (A; — j) 
avalanche following a j one. 

With the additional hypothesis of statistical independence between suc- 
cessive avalanches 

P{k-j\j)=p{k-j) (11) 

Defining C (x, /i) = (exp [ikx]) (the characteristic function) and Q (x, /i) = 
(7 (/c, /x) exp {ikx)) one obtains 

— C (x, /i) = 2Q (a;, /i) C (x, ji) — Q (x, /i) C (x, /i) (12) 

So far this is a very general equation that (with statistical independence of 
avalanches) would apply to any dynamical system. The dynamical informa- 
tion of each system is coded on the "merging kernel" '-f{k,fi). Using normal- 
ization (C (0,/i) = 1) and Kac's lemma (^—i-£:C (x,/i) \x=o= in Eq. lfT^ 
one obtains the following constraints on 'j{k,fj,) 



(7(fc,/i)) = 
{k-f{k,fi)) = 

A trivial solution to the constraints (fT^ is 



1 



(13) 



7(fc,/x) = - (14) 



leading to the equation 



fi^^P^ = C{x,f,)-C'{x,f,) (15) 

with solution 

Cix,f^) = —^ (16) 
H + q{x) 

with g (0) =0 and q (0) = —i, from normalization and Kac's lemma. Oth- 
erwise q (x) is an arbitrary function of x. 

However the Bak-Sneppen model is not in the class of solutions ()l(i|l . be- 
cause from this equation it follows that fx (^-citjt) ~ -'-) should be independent 
of /i, which is not verified by the numerical data (Fig. 5). 
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Figure 5: Real and imaginary parts of /i (^c^rj^ 



A power dependence of the merger kernel jlkjfi) 

7(/c, /i) = «! — h a2k^ 



together with the constraints ()T3|l implies 

which again, is not borne out by the data. 

However from the numerical simulation data {2ny = 2) for two close 
values of the measure and Eq. ()12p one may obtain direct information on 
7(/c, /i) from 

_ Jl^ Q (x, /i) exp {-ikx) dx 
f^^C {x, fi) exp {—ikx) dx 

An example is shown in Fig. 6, where both p{k,fi) and '-f{k,fi) are 
plotted. The general conclusion from these computations for several measures 
is that, starting from a non-zero value at = 1, 7 {k, fi) grows relatively fast 
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Figure 6: p{k) and ^{k) for n = 0.0022 



to the asymptotic value ^. A (weak) power dependence is only found for 
small k values. Therefore the behavior of the merger kernel, that codes the 
dynamical information of the Bak-Sneppen model, seems to be much more 
complex than thought before. 



4 Conclusions 

(i) The deterministic Bak-Sneppen model, studied in this paper, establishes 
a clear relation between the Lyapunov spectrum and scale-free behavior. 

Extremal behavior, where only a finite number of agents is active at each 
time, leads to the vanishing of all Lyapunov exponents in the infinite agents 
limit. On the other hand, local chaos insures that the convergence to zero is 
from above. These seem to be the most common ingredients in systems that 
display what has been called self-organized criticality. 

(ii) The self-organized state in the deterministic Bak-Sneppen model is 
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an interesting zero-measure subset of the invariant measure. Interpreting 
avalanches as return times to a set A, one concludes that a scaling law with 
exponent < 2 implies, by Kac's lemma, /i (A) = 0. Therefore one expects 
the zero measure feature to be present not only in the Bak-Sneppen model 
but also in other self-organized states with avalanche scaling laws. 

(iii) The avalanche equation, in the measure formulation discussed in 
Section 3, is a powerful guide to constraint the extraction of dynamical in- 
formation from dynamical systems of this type. In particular it provides pre- 
liminary evidence on a non-power behavior of the merger kernel. A stretched 
exponential behavior converging to - seems to be a better guess. This how- 
ever is a question that deserves further study. 



5 Appendix. Return times to small measure 



On a measurable space (M, B, /x), consider a differentiable map T for which 
H is an ergodic invariant measure. The return time ta {x) to a measurable 
set A starting from a point x e M is 



The probability to find a return time greater than A;, starting from an arbi- 
trary point X & M 



sets 



TA (x) ^ ini {{k > 1 \ T'^x e A) U oo 



(17) 




(18) 



and for a starting point x in A 



P {T^-^ > k) 



IJ'A {ta {x) > k) 



(19) 



where 





(20) 



Invariance of the measure implies 



ji {ta {x) > k) 



^x{A'^r^{TA{x)>k-l}) 

fj. {ta (x) > k-1)- II (A) Ha {ta (x) > k - 1) 



(21) 
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and because P (t'^^^ = k) = fi {ta (x) > k — 1) — fx {ta (x) > k) 

P (T^^^ = k)=^i{A)P (T^^^ > A; - 1) (22) 

Because of ergodicity, Poincare recurrence implies Xlfcli ^ [T^^^ = k^ = 1. 
On the other hand ^^=1 P {^^^^ > k - l) = kP (T^^^ = k) and for 

IX {A) 7^ one obtains Kac's lemma 

E {ta {x)\xeA)= / ta {x) dfiA (x) = (23) 

Defining, as in [H] 

c {k, A) = fiA {ta {x) >k)- 11 {ta (a;) > k) (24) 
one obtains from (PT|) and (fTT^ 

P (T^-^^ >k)={l-^i {A)) P (T^^^ > A: - 1) - Ac (A;, A) 
where Ac (fc, A) = c (A; — 1, A) — c (A;, A). By iteration 

P {T^-^ > A:) = (1 - (A))'^ - Ac (j, A){l-^^ {A)f-^ (25) 

i=i 

Let 

Ac (A) = sup I Ac (A;, A) I (26) 

k 

Then, by standard techniques, as in [Hj, one obtains the estimation 

|p^yA-.A ^ _g-MA)| < 2/i (A) + Ac (A) (1 - log Ac (A)) (27) 

This estimate is sharper than the one obtained in because it involves 
Ac (y4) rather than c{A) = supf. |c(A;,74)|. 

The estimate (jTTj) holds whenever /j, (A) ^ 0. However it is useless to find 
p{k, jj) = P (t^^^ = A;) in the /i [A) — > limit because 

_ ^~k^l{A) ^ ^-kt,{A) Cg-M(A) _i\ _^ 

Therefore, a different approach must be followed to control the yU [A) 
limit. 
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For general distributions of the form 

the characteristic function must satisfy the equation 

-i-C {x, fx) = i^C {x, fj,) + -C (x, //) 
dv ax II 

with solution 

f r dr dr \ 

C [x, II) = exp / — - - / — - 

Therefore a large family of solutions is obtained simply by specifying the 
function u (fi) (and its inverse ji (i/)). 
A particular case is 

Pa {k, n) ~ k-'^e-^'"^^!'^ 
For a — Q this is the geometric distribution with 

l^a {lA = -log(l - H) 

and for a = 1 

_ log(l -exp(ia;-t/)) 
log(l -exp(-i/)) 

with 

fi=il-eniog{l-e-') 
The functions Uq, (/x) for other values of a are plotted in Fig. 7. 
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Figure 7: Exponential argument ^'a(A*) for the probability densities 

Pa {k, jj) ~ exp {-kuo,{jJi)) 
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